#include "Voronoi.h"
#include <math.h>
#include <stdlib.h>
#include <iostream>
// Global Variables
vector<VoronoiRegion> vRegionList;
vector<DelaunayPoint> dPointList;

// Construct potential voronoi edges during voronoi region construction
// The line equation is in the form of y = ax + b
VoronoiEdge computeVoronoiEdge(DelaunayPoint p1, DelaunayPoint p2)
{
	double distance = sqrt(pow((p1.x-p2.x),2) + pow((p1.y-p2.y),2));
	double pii = pow(distance,2) - pow(p1.weight,2);
	double pij = pow(distance,2) - pow(p2.weight,2);

	bool slopeIsVertical;
	double a ;

	if(p1.y == p2.y){	slopeIsVertical = 1; a = 0;}
	else {				slopeIsVertical = 0; a = -(p2.x-p1.x)/(p2.y-p1.y);}

	Point midPoint;
	midPoint.x = p1.x - ((p1.x-p2.x) * pij/(pii+pij));
	midPoint.y = p1.y - ((p1.y-p2.y) * pij/(pii+pij));


	VoronoiEdge potentialEdge;
	potentialEdge.dp1 = p1;
	potentialEdge.dp2 = p2;
	potentialEdge.intersectionPoint = midPoint;
	potentialEdge.slope = a;
	potentialEdge.slopeIsVertical = slopeIsVertical;
	
	// two end points of voronoi edge
	VoronoiPoint vp1,vp2; // both infinity at the beginning
	
	if(!slopeIsVertical){
		vp1.x = midPoint.x-1, vp1.y = midPoint.y-a; // bottom end
		vp2.x = midPoint.x+1, vp2.y = midPoint.y+a;} // upper end
	else{
		vp1.x = midPoint.x, vp1.y = midPoint.y-1; // bottom end
		vp2.x = midPoint.x, vp2.y = midPoint.y+1;
	}

	vp1.isInfinity = true;
	vp2.isInfinity = true;

	// assign corresponding delaunay edge
	vp1.dp1 = p1, vp1.dp2 = p2;
	vp2.dp1 = p1, vp2.dp2 = p2;
	// vp1.degree = 2, vp2.degree = 2;

	potentialEdge.vp1 = vp1;
	potentialEdge.vp2 = vp2;
	 
	return potentialEdge; // will be destroyed or not ?
}

Point calcIntersection(VoronoiEdge ve1, VoronoiEdge ve2, bool& isIntersect)
{
	Point p;
	//ve1C is the intersection at y-axis
	if(!ve1.slopeIsVertical && !ve2.slopeIsVertical){
		double ve1C = ve1.intersectionPoint.y - ve1.slope*ve1.intersectionPoint.x;// c = y -mx
		double ve2C = ve2.intersectionPoint.y - ve2.slope*ve2.intersectionPoint.x;
		p.x = (ve2C-ve1C) / (ve1.slope-ve2.slope);
		p.y = ve1.slope*p.x + ve1C;
	}
	else if(ve1.slopeIsVertical && !ve2.slopeIsVertical){
		double ve2C = ve2.intersectionPoint.y - ve2.slope*ve2.intersectionPoint.x;
		p.x = ve1.intersectionPoint.x;
		p.y = ve2.slope*p.x + ve2C;
	}
	else if(ve2.slopeIsVertical && !ve1.slopeIsVertical){
		double ve1C = ve1.intersectionPoint.y - ve1.slope*ve1.intersectionPoint.x;
		p.x = ve2.intersectionPoint.x;
		p.y = ve1.slope*p.x + ve1C;
	}
	else if(ve1.slopeIsVertical && ve2.slopeIsVertical){
		isIntersect = 0;
	}
	
	int ve1intersected,ve2intersected;
	
	if(ve1.slope == ve2.slope && ve1.slopeIsVertical == 0 && ve2.slopeIsVertical == 0){
		isIntersect = 0;
	}
	else{
		if		(ve1.vp1.isInfinity && ve1.vp2.isInfinity)ve1intersected = 1;
		else if	(ve1.vp1.isInfinity && !ve1.vp2.isInfinity){
			if(ve1.vp2.x > ve1.vp1.x){
				if(p.x < ve1.vp2.x) ve1intersected = 1;
				else ve1intersected = 0;
			}
			else if (ve1.vp2.x < ve1.vp1.x){
				if(p.x > ve1.vp2.x) ve1intersected = 1;
				else ve1intersected = 0;
			}
		}
		else if	(!ve1.vp1.isInfinity && ve1.vp2.isInfinity){
			if(ve1.vp1.x > ve1.vp2.x){
				if(p.x < ve1.vp1.x) ve1intersected = 1;
				else ve1intersected = 0;
			}
			else if (ve1.vp1.x < ve1.vp2.x){
				if(p.x > ve1.vp1.x) ve1intersected = 1;
				else ve1intersected = 0;
			}
		}
		else if	(!ve1.vp1.isInfinity && !ve1.vp2.isInfinity){
			if((ve1.vp1.x>p.x && p.x>ve1.vp2.x) || (ve1.vp2.x>p.x && p.x>ve1.vp1.x))ve1intersected = 1;
			else ve1intersected = 0;
		}
		
		if(ve2.vp1.isInfinity && ve2.vp2.isInfinity)ve2intersected = 1;
		else if	(ve2.vp1.isInfinity && !ve2.vp2.isInfinity){
			if(ve2.vp2.x > ve2.vp1.x){
				if(p.x < ve2.vp2.x) ve2intersected = 1;
				else ve2intersected = 0;
			}
			else if (ve2.vp2.x < ve2.vp1.x){
				if(p.x > ve2.vp2.x) ve2intersected = 1;
				else ve2intersected = 0;
			}
		}
		else if	(!ve2.vp1.isInfinity && ve2.vp2.isInfinity){
			if(ve2.vp1.x > ve2.vp2.x){
				if(p.x < ve2.vp1.x) ve2intersected = 1;
				else ve2intersected = 0;
			}
			else if (ve2.vp1.x < ve2.vp2.x){
				if(p.x > ve2.vp1.x) ve2intersected = 1;
				else ve2intersected = 0;
			}
		}
		else if(!ve2.vp1.isInfinity && !ve2.vp2.isInfinity){
			if((ve2.vp1.x>p.x && p.x>ve2.vp2.x) || (ve2.vp2.x>p.x && p.x>ve2.vp1.x))ve2intersected = 1;
			else ve2intersected = 0;
		}

		if(!ve1intersected || !ve2intersected){
			isIntersect = 0;
		}
	}

	return p;
}

double dotProduct(Point a, Point b)
{
	return a.x*b.x + a.y*b.y;
}

void updateVoronoiEdge(VoronoiEdge& ve1, VoronoiEdge& ve2, DelaunayPoint cp, Point p)
{
	Point vector1, vector2, vector3, vector4, vector5;
	vector1.x = cp.x - p.x;
	vector1.y = cp.y - p.y;

	vector2.x = ve1.vp1.x - p.x;
	vector2.y = ve1.vp1.y - p.y;
	vector3.x = ve1.vp2.x - p.x;
	vector3.y = ve1.vp2.y - p.y;

	vector4.x = ve2.vp1.x - p.x;
	vector4.y = ve2.vp1.y - p.y;
	vector5.x = ve2.vp2.x - p.x;
	vector5.y = ve2.vp2.y - p.y;

	if ( dotProduct(vector2, vector4) > 0 )
	{
		if ( dotProduct(vector2, vector1) >
			 dotProduct(vector4, vector1) )
		{
			ve1.vp2.x = p.x;
			ve1.vp2.y = p.y;
			// to add ve1.vp2.dp3
			ve1.vp2.isInfinity = false;

			ve2.vp1.x = p.x;
			ve2.vp1.y = p.y;
			// to add ve2.vp1.dp3
			ve2.vp1.isInfinity = false;
		}
		else
		{
			ve1.vp1.x = p.x;
			ve1.vp1.y = p.y;
			// to add ve1.vp2.dp3
			ve1.vp1.isInfinity = false;

			ve2.vp2.x = p.x;
			ve2.vp2.y = p.y;
			// to add ve2.vp1.dp3
			ve2.vp2.isInfinity = false;
		}
	}
	else
	{
		if ( dotProduct(vector2, vector1) >
			 dotProduct(vector5, vector1) )
		{
			ve1.vp2.x = p.x;
			ve1.vp2.y = p.y;
			// to add ve1.vp2.dp3
			ve1.vp2.isInfinity = false;

			ve2.vp2.x = p.x;
			ve2.vp2.y = p.y;
			// to add ve2.vp1.dp3
			ve2.vp2.isInfinity = false;
		}
		else
		{
			ve1.vp1.x = p.x;
			ve1.vp1.y = p.y;
			// to add ve1.vp2.dp3
			ve1.vp1.isInfinity = false;

			ve2.vp1.x = p.x;
			ve2.vp1.y = p.y;
			// to add ve2.vp1.dp3
			ve2.vp2.isInfinity = false;
		}
	}
}

bool checkDeleteEdge(VoronoiEdge target, VoronoiEdge newEdge, DelaunayPoint cp)
{
	Point vector1, vector2, vector3;
	vector1.x = cp.x - newEdge.intersectionPoint.x;
	vector1.y = cp.y - newEdge.intersectionPoint.y;


	// parallel case hasn't been considered
	vector2.x = target.vp1.x - newEdge.intersectionPoint.x;
	vector2.y = target.vp1.y - newEdge.intersectionPoint.y;
	vector3.x = target.vp2.x - newEdge.intersectionPoint.x;
	vector3.y = target.vp2.y - newEdge.intersectionPoint.y;

	if (dotProduct(vector1, vector2)>0 && dotProduct(vector1, vector3)>0)
		return false;
	else
		return true;

}

void constructVoronoiRegion(DelaunayPoint cp, VoronoiRegion& vRegion)
{
	vRegion.dp1.x = cp.x;
	vRegion.dp1.y = cp.y;

	vector<DelaunayPoint>::iterator it = dPointList.begin();
	while ( it != dPointList.end() )
	{
		if ( it->x != cp.x || it->y != cp.y )
		{
			VoronoiEdge newEdge = computeVoronoiEdge(cp, *it);
			vector<VoronoiEdge>::iterator it2 = vRegion.vEdgeList.begin();
			int count = 0;
			if ( vRegion.vEdgeList.size() == 0 )
			{
				vRegion.vEdgeList.push_back(newEdge);
				continue;
			}		

			while ( it2 != vRegion.vEdgeList.end() )
			{
				bool isIntersect;
				Point p = calcIntersection(newEdge, *it2, isIntersect);
				
				if ( isIntersect )
				{
					updateVoronoiEdge(newEdge, *it2, cp, p);
					count++;
				}
				else
					checkDeleteEdge(*it2, newEdge, cp);	
				it2++;
			}
			if (count != 0)
				vRegion.vEdgeList.push_back(newEdge);
		}
		it++;
	}

}

void constructVoronoiDiagram(void)
{
	vector<DelaunayPoint>::iterator it = dPointList.begin();
	while ( it != dPointList.end() )
	{ 
		VoronoiRegion newVRegion;
		newVRegion.dp1.x = it->x;
		newVRegion.dp1.y = it->y;
		constructVoronoiRegion(*it, newVRegion);
		it++;
	}
}

int main(){
	DelaunayPoint dp1,dp2,dp3;

	dp1.x = 0;
	dp1.y = 0;
	dp1.weight = 2;
	
	dp2.x =5;
	dp2.y =0;
	dp2.weight = 2;

	dp3.x =0;
	dp3.y =-4;
	dp3.weight = 2;

	VoronoiEdge ve1 = computeVoronoiEdge(dp1,dp2);
	VoronoiEdge ve2 = computeVoronoiEdge(dp2,dp3);
	VoronoiEdge ve3 = computeVoronoiEdge(dp3,dp1);

	cout<< dp1.x <<  " " << dp1.y << "\n";
	cout<< dp2.x <<  " " << dp2.y << "\n";
	cout<< dp3.x <<  " " << dp3.y << "\n";

	cout<< ve1.intersectionPoint.x << " " << ve1.intersectionPoint.y <<  " " << ve1.slope << " " << ve1.vp1.x << " " << ve1.vp1.y << " " << ve1.vp2.x << " " << ve1.vp2.y << "\n";
	cout<< ve2.intersectionPoint.x << " " << ve2.intersectionPoint.y <<  " " << ve2.slope << " " << ve2.vp1.x << " " << ve2.vp1.y << " " << ve2.vp2.x << " " << ve2.vp2.y << "\n";
	cout<< ve3.intersectionPoint.x << " " << ve3.intersectionPoint.y <<  " " << ve3.slope << " " << ve3.vp1.x << " " << ve3.vp1.y << " " << ve3.vp2.x << " " << ve3.vp2.y << "\n";
	
	//bool a,b,c;
	//Point int1 = calcIntersection(ve1,ve2,a);
	//Point int2 = calcIntersection(ve2,ve3,b);
	//Point int3 = calcIntersection(ve3,ve1,c);

	return(0);



}